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Introduction 

The  accuracy  of  the  numerical  solution  of  a Partial  Differential  Equation  (PDE)  is  a function  of  the 
order  of  the  leading  terms  of  the  truncation  error  of  the  numerical  algorithm.  These  terms  in  turn 
depend  both  on  the  higher  order  derivatives  of  the  solution  variables  and  on  the  local  mesh  size.  In 
practical  applications  it  is  important  to  guarantee  a sufficient  level  of  accuracy,  if  not  mesh 
independence,  of  the  solution.  This  condition  is  normally  reached  by  comparing  solutions  obtained 
on  successively  refined  grids,  a procedure  that  is  empirical  and  time-consuming,  and  further 
complicated  by  the  fact  that  mesh  generation  is  conducted  heuristically,  since  the  solution  is  not 
known  in  advance. 

Given  the  nature  of  the  truncation  error  of  the  numerical  algorithms  for  the  solution  of  PDEs,  two 
basic  methods  for  improving  the  accuracy  of  a solution  can  be  readily  identified:  a)  construction  of 
higher-order  schemes  or,  given  a fixed  numerical  scheme,  b)  to  control  the  numerical  error  via  grid 
optimization. 

In  industry,  while  run-time  efficiency  is  important,  access  to  a code’s  source  may  not  be  guaranteed 
and  it  can  be  impractical  to  modify  existing  codes.  In  addition,  commercial  code  vendors  may  not 
be  over-enthusiastic  about  introducing  a generic  and  revolutionary  mesh  adaptation  process, 
because  it  may  show  the  shortcomings  of  the  current  version  of  their  code.  Therefore,  a more 
effective  approach  consists  in  adapting  the  grid  to  the  solution,  without  access  to  the  source.  The 
process  ideally  begins  with  an  arbitrary  mesh,  hopefully  coarse  and  generated  by  a non-expert,  and 
leads  to  a highly  adapted  customized  grid,  tailor-made  to  each  application.  Thus,  no  expensive 
algorithm  modifications  are  required  and  any  CFD  code  can  satisfy  a predefined  level  of  accuracy, 
as  desired  or  as  required  by  the  problem  specifications,  without  user  intervention.  This  strategy  is 
ideal  for  parameter  sensitivity  studies,  for  example,  where  all  the  different  solutions  could  be 
obtained  on  automatically  adapted  grids  generated  from  a single  original  grid.  In  such  case  one 
becomes  sure  that  the  change  in  performance  I sure  to  the  change  in  parameter  and  not  due  to  any 
sensitivity  to  the  grid. 

The  code  discussed  here  represents  an  advanced  novel  anisotropic  mesh  adaptation  approach  for 
unstructured  grids  that  allows  the  customization  of  grids  based  on  a specified  error  level.  Equally 
important,  it  maintains  the  CAD  integrity  of  the  given  geometry.  Perhaps  the  best  showcase  of  the 
abilities  of  this  software  is  the  adaptation  of  unstructured  tetrahedral  meshes  over  a complete 
Airbus  320  and  a Boeing  747  with  flow-through  engine  nacelles. 


Paper  presented  at  the  RTO  AVT  Symposium  on  “ Reduction  of  Military  Vehicle 
Acquisition  Time  and  Cost  through  Advanced  Modelling  and  Virtual  Simulation”, 
held  in  Paris,  France,  22-25  April  2002,  and  published  in  RTO-MP-089. 
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Methodology 

Any  (anisotropic)  mesh  adaptation  algorithm  requires  two  fundamental  components:  an  a- 
posteriori  error  estimator  and  a mechanism  to  coarsen,  refine,  swap  and  stretch  the  mesh.  In 
general  it  is  impossible  to  evaluate  the  exact  truncation  error  of  the  numerical  approximation  of  a 
PDE  or  system  of  PDEs,  since  the  numerical  approximation  is  only  piecewise  continuous  and  in  any 
case  the  exact  solution  is  unknown.  Most  of  the  error  estimators  in  use  are  based  on  the  first 
derivative  or  the  residual  of  a user-specified  solution  variable.  The  error  estimator  proposed  here  is 
based  on  the  Hessian  (matrix  of  second  derivatives),  which  is  a better  representation  of  the 
truncation  error  between  the  PDEs  and  the  discrete  equations  actually  solved.  The  Hessian  matrix 
provides  not  only  information  on  the  error  density,  but  also  on  its  directionality,  thus  making 
anisotropic  mesh  adaptation  possible. 


The  edge-based  error  density  yean  be  written  as: 


7= 


As  M As  dt 


where  As  is  the  edge  length,  M is  the  Hessian  matrix  and  t is  the  non-dimensional  parametrization 
variable  along  the  edge.  Once  the  error  density  is  determined  for  every  edge  in  the  grid,  the  user 
must  select  the  desired  error  level  for  the  adaptation.  Unless  dictated  by  the  specific  problem,  an 
appropriate  error  density  can  be  selected  visually  with  any  graphic  interface. 

The  Mesh  adaptation  procedure  consists  of  four  main  steps:  a)  edge  coarsening,  b)  edge  refinement, 
c)  edge  swapping  and  d)  node  movement.  Steps  a)  and  b)  can  be  deemed  binary  operations,  since 
two  connected  edges  are  combined  into  one,  or  one  edge  is  split  in  two,  respectively.  These  two 
operations  may  produce  grids  that  are  not  optimal,  therefore  the  edge  swapping  mechanism  will 
correct  the  orientation  of  any  edge  that  might  not  be  properly  positioned  and  lead  to  distorted 
meshes.  Node  movement  stretches  and  realigns  the  edges  along  the  proper  direction,  and  ensures 
that  the  error  is  equi-distributed  and  the  maximum  and  minimum  edge  lengths  and  the  minimum 
aspect  ratio  are  respected.  Coarsening,  refinement  and  swapping  can  be  regarded  as  accelerators 
for  mesh  movement. 


Node  movement  is  based  on  a spring  analogy,  where  each  element  edge  in  the  grid  is  represented  by 
a spring  with  a coefficient  of  elasticity  proportional  to  the  error  density.  The  goal  in  this  case  is  to 
equi-distribute  the  error  across  the  computational  domain.  The  original  CAD  surface  definitions 
are  used  at  each  step  to  ensure  that  surface  nodes  are  reprojected  on  the  proper  surfaces  and  the 
original  geometry  is  fully  respected.  Since  the  mesh  adaptation  procedure  is  a nonlinear  process, 
steps  a)  to  d)  are  repeated  in  cycles  at  every  adaptation. 

Several  control  parameters  are  available  to  control  the  quality  of  the  final  mesh,  among  them 
aspect  ratio  and  a maximum  deviation  on  curvature.  High  maximum  aspect  ratio  values  can  lead  to 
extremely  stretched  grids,  and  while  FEN SAP  has  no  difficulties  handling  highly  stretched  grids, 
other  flow  solvers  might.  Maximum  deviation  on  curvature  allows  one  to  determine  the  acceptable 
level  of  coarseness  on  curved  surfaces  that  prevents  poor  surface  rendition.  This  may  be  noticeable 
on  the  lips  of  the  nacelles  or  on  wing  tips,  where  surface  curvature  is  high.  In  general  the  meshes 
adapted  show  much  more  faithful  renditions  of  the  original  CAD  than  even  the  finest  original 
meshes. 


There  is  no  question  that  mesh  adaptation  provides  a useful  and  necessary  tool.  The  expense 
incurred  in  adapting  the  mesh  and  recomputing  the  flow  solution  again  and  again  is  not  negligible, 
therefore  the  focus  of  our  research  is  to  provide  some  guidelines  to  facilitate  the  choice  of 
parameters  by  illustrating  their  effect  on  the  mesh  quality  and  overall  solution  time.  The  subject  of 
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grid  convergence,  important  from  the  point  of  view  of  computational  efficiency,  will  also  be 
addressed. 

Results:  Mesh  Adaptation  on  a Boeing  747  and  an  Airbus  320 

Both  the  Boeing  747  and  Airbus  320  feature  fuselage-wing-pylon-nacelle  geometries,  but  the  Airbus 
320  has  simple  flow-through  engine  nacelles,  while  the  Boeing  747  has  both  engine  nacelles  and  fan 
shrouds  and  in  addition  some  detail  of  the  cockpit  windows.  For  both  geometries,  the  outer 
boundary  is  a hemisphere  with  a diameter  of  approximately  40  fuselage  lengths.  Flight  conditions 
for  the  Boeing  747  are:  Mach  number=0.88,  altitude  of  35,000ft  and  AoA  of  4°.  Flight  conditions  for 
the  Airbus  are:  Mach  number=0.8,  altitude  of  30,000ft  and  AoA  of  3°.  The  CFD  solutions  were 
obtained  with  FENSAP  (Finite  Element  Navier-Stokes  Analysis  Package),  also  a proprietary 
technology  of  Numerical  Technologies  International. 

For  the  Boeing  747,  the  initial  tetrahedral  grid  contains  431,181  nodes  and  2.3  million  linear 
tetrahedral  elements.  The  error  density  target  value  was  set  y=0.15  and  remained  constant  for  all 
the  subsequent  adaptations.  The  minimum  and  maximum  element  edge  sizes  were  set  at  0.025  in 
and  25,000  inches,  respectively,  compared  to  a fuselage  length  of  approximately  2,700  inches.  Eight 
complete  cycles  of  solution  and  adaptation  were  performed.  After  six  adaptations,  the  final  grid  size 
had  decreased  to  242,310  nodes  and  1.3  million  tetrahedral  elements,  an  decrease  of  43%,  for  an 
immense  improvement  in  accuracy. 

It  should  be  noted  that  some  intermediate  grids,  say  from  iteration  1-3,  may  become  somewhat 
larger  than  the  original  as  mesh  adaptation  begins  to  cluster  the  mesh  near  the  salient  features 
(discontinuities  and  vortices  in  this  case).  This  phenomenon  is  related  to  the  density  of  the  original 
mesh,  and  how  well  it  captures  the  initial  solution.  Initially  the  discontinuities  are  smeared  and 
require  much  refinement,  but  as  they  become  sharper  stretching  and  coarsening  remove  some  of 
the  grid  points  and  the  mesh  size  begins  to  decrease.  Note  that  in  3D  discontinuities  are  surfaces 
and  hence  it  is  difficult  to  estimate  the  extent  of  the  refinement  just  by  looking  at  the  surface  mesh. 

The  adapted  grid  shows  substantial  coarsening  on  the  lower  wing  surface  and  refinement  on  the 
upper  surface.  The  outboard  wing  shock  is  very  well  captured,  and  the  inboard  shock  is  also 
beginning  to  appear.  A comparison  of  the  initial  and  final  upper  wing  grids  and  corresponding 
solutions  of  the  Airbus  320  are  shown  in  Figure  1.  Figures  2,  3 and  4 show  some  surface  grid  details 
and  solutions  of  the  Boeing  747.  Note  in  Figure  5 how  adaptation  has  improved  the  discretization  of 
the  fan  cowl  lip.  Figure  6 shows  the  Cp  distribution  at  three  spanwise  locations  before  and  after 
adaptation.  Finally,  Figure  7 shows  the  sequence  of  error  distributions  in  the  computational  domain 
as  the  adaptation-solution  cycles  evolve.  Note  that  there  is  little  change  from  cycle  5 to  cycle  6 and 
the  adaptation  procedure  is  converged. 


Conclusions 

Mesh  adaptation  is  a powerful  and  cost-effective  tool  that  can  be  used  to  increase  (or  decrease  if  so 
desired),  to  a specified  and  uniform  level,  the  computational  accuracy  of  any  given  CFD  code 
without  extensive  algorithm  redesign  or  initial  investment  in  grid  generation.  Furthermore  the  cost 
of  ensuring  mesh  independence  is  substantially  reduced,  and  the  error  level  of  the  results  can 
always  be  guaranteed. 
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Figure  1.  Airbus  320:  initial  and  final  grids  and  solutions 
after  6 adaptation  cycles 


Figure  2.  Boeing  747:  initial  and  final  grids  after  6 adaptation  cycles 
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Figure  3.  Boeing  747:  initial  and  final  grid,  pylon-nacelle  detail 


Figure  4.  Boeing  747:  initial  and  final  grid,  pylon-nacelle  detail 


Figure  5.  Boeing  747:  initial  and  final  grid,  fan  cowl  lip  detail 
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Figure  6.  Boeing  747:  initial  and  final  Cp  distribution  at  three  wingspan  stations 


